Electrostatic interactions play an essential role in molecular and cellular processes that include signal transmission at synaptic junctions, ion-transport, molecular recognition, and stability and function of DNA, RNA and proteins. Of paramount biological importance is the collective behavior of ions, molecules, and macromolecules having inhomogeneous charge distributions in the aqueous crowded environment of a living cell. Central to this environment is water, which is a complex solvent with non-bulk properties near ions, molecules and water-biopolymer interfaces. The main challenge for modeling electrostatics, as the core calculation of Molecular Dynamics and Monte Carlo simulations, is to balance the accuracy of interactions among charges and the efficiency required for realistic biological systems. For all-atom simulations the Ewald method is currently the most accurate method for large enough system where periodic boundary effects are negligible. Unfortunately, the simulation of a macromolecule requires a large simulation volume, making the Ewald method and similar methods too expensive while a smaller system size compromises accuracy. Meanwhile, modeling the solvent using continuum electrostatics is much faster, but accuracy is compromised by neglecting molecular scale inhomogeneity, such as atomic details near the surface of the solute.